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Abstract 



We study the stationary state of a simple exclusion process on a ring which was 
recently introduced by Arndt et al [ll^l^]. This model exhibits spatial condensation 



of particles. It has been argued |14, 15] that the model has a phase transition from 
a "mixed phase" to a "disordered phase". However, in this paper exact calculations 
are presented which, we believe, show that in the framework of a grand canonical 
ensemble there is no such phase transition. An analysis of the fluctuations in the 
particle density strongly suggests that the same result also holds for the canonical 
ensemble. 



1 Introduction 

One-dimensional driven diffusive systems have attracted the interest of researchers in a 
wide variety of fields, including mathematics, physics, chemistry and biology . One of 
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the most basic of these models is the asymmetric simple exclusion process (ASEP), which 
already exhibits a rich behavior. Within physics, for example, the ASEP has been used to 
study boundary induced phase transition [^, spontaneous symmetry breaking 0, phase 
separation and shocks . As one of the few models in non-equilibrium physics which 
is analytically accessible, it has served as a laboratory to test basic concepts of statistical 
mechanics for far from equilibrium systems. Moreover, its simplicity has led to applications 
in such diverse fields as the kinetics of biopolymerization fl^ , vehicular traffic |[ll|,|l^ and, 
most recently, sequence alignment problems in biology [0. 

In this paper, we study a related exclusion process which exhibits a spatial condensation 
of particles. This model was introduced by Arndt, Heinzel and Rittenberg [|14|,|15| and will 
be called AHR model below. 

In the AHR model two types of particles, called positive and negative, occupy the sites 
of a periodic one- dimensional lattice of length L. The particles are subject to hard-core 
exclusion, so that there are three possible states at each site: empty, occupied by a positive 
particle, or occupied by a negative particle. Positive particles hop at rate a into empty 
spaces (holes) to their immediate right, and negative particles at the same rate into empty 
spaces to their left. Adjacent positive and negative particles also exchange positions with 
asymmetric rates pr and pl- Specifically, the rules are 

+ + with rate a, 

— — with rate a, 

H — — h with rate pr, 

— h — > H — with rate p^, 

we assume that pr, a > and > 0. A key parameter in the model is ratio of pl to pr, 



which we denote by g: q = Pl/pr- Our notation here is not that of |jTj,|T5[, where q denotes 
the inverse of this ratio; we adopt this change to agree with the standard notation in the 
theory of basic hypergeometric functions (see e.g. Finally, we normalize our time 

scale by setting pr = 1, so that the rates for exchange of positive and negative particles 
are 1 and q. Note that the numbers N± of positive and negative particles are constants of 
the dynamics; we will write p± = N±/L for the particle densities. 

We remark that when Pl = and pr = a = 1 the model corresponds to a two-species 
ASEP (with "first" and "second class" particles) which has been solved in |^. The system 
with no holes (A^+ + = L) corresponds to the partially asymmetric exclusion process 
on a ring, for which the stationary measure is trivial (all configurations are equally likely). 

The stationary state of the AHR model displays a rich dynamical behavior, which was 



extensively investigated in Two methods were used: Monte Carlo simulations were 

carried out in the canonical ensemble (CE), in which the numbers N± of particles are fixed, 
and analytic calculations (still in finite systems) using a quadratic algebra were carried out 
in a grand canonical ensemble (GCE), in which the densities are controlled by fugacities 
^± and can fiuctuate. These results were then extrapolated to infinite volume. Most of the 
work was done in the case of equal densities {N^ = or = C-); which is the only case 
we will consider here. 
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These investigations strongly suggest the existence of three different phases for the 
model, differing in the nature of typical spatial particle configurations (see Figure 1 for 
space time representations). We summarize the discussion of [H, 15| and the nomencla- 
ture adopted there. For g > 1 the system is in the "pure phase": typical configurations 
contain three blocks, one each of holes, positive particles, and negative particles. In the 
thermodynamic limit (L —>■ oo) the current vanishes. The system in this phase is closely 
related to a class of models studied in 0. For < g < 1 the phase depends on the density 
p = p_i_ = p_. For relatively small densities, density profiles of both species are uniform 
and there is no spatial condensation (the "disordered phase"). For larger densities one 
observes the formation of a condensate made of + and — particles (with some holes at the 
boundaries); a block of holes, occasionally traversed by + or — particles, occupies the rest 
of the system (the "mixed phase" ) . At the transition between the mixed and disordered 
phases the Monte Carlo simulations suggest that the infinite volume current J of + parti- 
cles increases monotonically from zero to a value (1 — g)/4 as the the density p increases 
from zero to a certain value pc, and then takes the constant value (1 — g)/4 for < P < 1/2 
(see Figure 8 of [|I^). This indicates that the function J(p) is non-analytic — in fact, the 
simulations suggest it is not differentiable — at a certain value of density p^. The same 
apparent non-analyticity is seen in plots, derived from calculations in the GCE, of J as 
a function of q at fixed p (Figures 1 and 2 of [0). This would correspond to a "phase 
transition" in the language of equilibrium physics. 

In fact, however, we believe that our analysis shows that these numerical results are 
misleading. Our main result is an explicit exact formula for the current J in the infinite 
system in the framework of a GCE. This computation was possible because the AHR model 
is closely related to the one- dimensional partially asymmetric simple exclusion processes, 
for which recently analytic exact expressions for various expectation values have been 
computed ||T7Hl9[. It turns out that if g < 1 then the current J(p) is analytic everywhere, 
although the derivative of J(p) changes very rapidly as p passes pc when q is close to one. 
This strongly suggests that in the GCE, there is no phase transition. For this reason, we 
will here speak of the mixed regime and the disordered regime, and refer to the union of 
these regimes as the mixed/disordered phase. 

Now, is this also true for the CE? Based on the observation that the density fluctuations 
in the GCE go to zero in the thermodynamic limit, we believe that the answer to this 
question is yes, although the argument is not rigorous (See section 5 for a more detailed 
discussion.) 

The existence of a phase transition between the pure and the mixed/disordered phase 
seems uncontroversial |1^,|15| and will not be discussed in the following. Therefore, q is 
normally chosen to satisfy < g < 1 from now on, although we will occasionally consider 
g = 1 as a special limiting case. 

This paper is organized as follows. In section 2, we present the basic structure of 
the exact solution and define the GCE and CE. The current and density in the GCE are 
computed in section 3. An analysis of the model when g is close to one is carried out in 
section 4. In section 5 we discuss to what extent the results for the GCE are valid for the 
CE. The final section contains our conclusions. 
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2 Exact solution of the stationary state 



The stationary state of the AHR model can be constructed by applying the so-called matrix 
product ansatz (MPA) |^,^; this is the method used for the grand canonical calculations 





various exact results for 
and references 
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in IjTj, 0. Since its original application to the ASEP in 
exclusion processes have been achieved via the MPA; see for example 
therein. Since the method is now well known, we will skip a detailed presentation of the 
formalism and even omit the proof that the particular version presented here correctly 
determines the stationary state. 

The MPA expresses the probability of any configuration in the stationary state of the 
AHR model, up to an overall normalization, as a trace of a product of certain matrices. 
Thus, for example, the probability P of the configuration H — 00+ can be written in the 
form P(H — 00+) = Z^^Tt{DEAAD); here Z is an ensemble-dependent normalization 
constant which will be discussed shortly. The matrices D, E and A are assigned to the 
local states +, — and 0, respectively, and are required to satisfy the following algebraic 
conditions, which guarantee stationarity of the resulting state: 



DE-qED = C{D + E), 
aDA = (A, 
aAE = (A. 



(2.1) 



Here the number ( is arbitrary; a natural choice, which we adopt from now on, is 

C = 1 - g. (2.2) 
The conditions (|2.1| ) are closely related to those for the matrices of the MPA for the 



partially asymmetric simple exclusion processes |T7HT9[] 
It is convenient to introduce the new parameter 



-1 + 



1-q 



a 



(2.3) 



In terms of a, explicit representations for the matrices D, E and A which satisfy (|2.1|) are 
given by 



D 



1 + a 





/Ci 

1 + aq 






1 + 








and 



E = D^, 
A=\0){0\, 



(2.4) 
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where the superscript T indicates the transpose, 



c„ = (l-g")(l-aV"^), (2.5) 

(0| = [1 ... ], and |0) = (Of. 

Recall that we take < g < 1, and certainly a > —1. Thus if a < 1 then c„ is 
positive for any n (> 1), and the roots in (|2.4|) are real. If a > 1 then the roots may 
be imaginary; moreover, c„ can vanish for special values of a and n, leading to finite 
dimensional representations. The case a > 1 requires special treatment in other ways, 
which we will mention briefly in section 3. In the future, unless special mention is made, 
we assume that — 1 < a < 1; this is the interesting case since the transition between the 



mixed and disordered regime occurs in this region (see Table 1 of [|15|; our results below 
show that this transition in fact occurs for a < 0). 

Since the matrices in (|2.4|) are infinite dimensional, it is not a priori clear that the 
trace of an arbitrary product of these matrices exists. However, A is a projector of rank 
one, so that the trace of any product of these matrices will exist as long as there is at 
least one A matrix in the product or, equivalently, as long as there is at least one hole in 
the corresponding configuration. We will consider here only such well defined traces since, 
as mentioned in the introduction, when there are no holes the model has a trivial steady 
state. 

For the CE, in which are specified, the normalization constant Z2\^ is the sum 
of traces of matrix products over all configurations which respect the given A^_|_ and A^_. 
This ensemble is in many ways the natural one; Monte Carlo simulations for example, are 
usually done in the CE. However, it seems difficult to compute It is easier to 

analyze the model in a superposition of CE's, that is, in an appropriately defined GCE. 
There is a certain ambiguity in the definition of such an ensemble, but a standard method 



has emerged P,[T5|. In the AHR model it would be natural to introduce different fugacities 
^± for the two types of particles, but since we are interested in the case of equal densities, 
we set ^+ = ^_ = ^. Thus we take the probability of a configuration T1T2 ■ ■ ■ tl, where 
Tj = +, — or for j = 1, 2, . . . , L, to be 

1 ^ 

P«c^(n, r2, . . . , n) = ^Tr l[[6r„+^D + d^^^^E + 5.„o^]- (2.6) 



n=l 



Now the normalization constant Zl{^) = Z^^^(^) is a sum over traces of matrix products 
for all configurations of size L which have at least one hole, that is, 

Zl{0 = Tr'G^ (2.7) 

with 

G = A + ^C (2.8) 

where 

C = D + E. (2.9) 
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Here the prime in Tr' indicates a sum over all terms with at least one A. A more detailed 
discussion of the relationship between the CE and the GCE will be given in section 5. 

The current 7+ of + particles in the GCE is defined as the average flux of + particles 
through a given bond. It is the same as the current J_ of — particles, defined similarly. 
The current J — J+ — J- and the density p = p+ = p- are easily seen to be given by 

MO = (1 - (2-10) 

PLiO = (2.11) 

where the subscript L indicates the system size. As pointed out in the introduction, wc are 
interested in the infinite volume limits J((^) = limL^oo Jl{0 and p{0 = limL^ooPL({) of 
these quantities because a phase transition between the mixed and the disordered regime 
should manifest itself as non-analyticity in J as a function of p. J(^) and p(^) will be 
computed in the next section. 



3 Exact Results in the Thermodynamic Limit 

Instead of considering directly, we introduce the generating function 

oo 

e(A,o = $]A^-^z40- (3.1) 

L=l 

This sum can be explicitly evaluated as follows. First we rewrite ^l(0 as a summation of 
the terms with the condition that the left- most hole is at site j + 1: 

L-1 L-1 

Zl{0 = J2TT{^CyAG'^-^-' = 5](0|G^-^-i(eC)^|0), (3.2) 

j=0 j=0 

where we used the cyclicity of trace and the explicit representation of = |0)(0|. Hence 
the generating function 0(A, ^) is rewritten as 

oo 

0(A,O= E(0|(^^)'(^^^)'|0)- (3-3) 

j,k=0 

Second we expand as 

G'^^Yl Yl iccy°A'^'{^cy' ■■■A"'^{ccy^. (3.4) 

r>0 . . . „ 

:;o,>>o 

ji,mi,... ,jr-l,mr>0 
jo+mi+jiH [-mr+jr=k 
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(3.5) 



If we observe that = A and define 

oo 
n=0 

so that 

Y.m'{o\no) = (^) , Y.{n + mxr (omo) = -^^x (^) , (3.6) 

then after some computation we obtain 

e(A,o 



(3.7) 



Finally it turns out that the function x{^) has been known in mathematics litera- 
ture [^. This is related to the fact that the matrix C is essentially the Jacobi matrix as- 



sociated with certain g-orthogonal polynomials, called Al-Salam-Chihara polynomials [23 



(and regarded as a special case of the Askey- Wilson polynomials p4|); see [0 for an 
explanation of this connection. For example, from [Q, equation (3.21), it follows that 



with 



xix) = f{y{x)) 



r( . (?z/^g)oo(g;g)oo 



(ay; q) 



loo 



y{x) 



{ay, q) 

X — 2 — — 4x 



n=0 



{qy^;<i)n{q; q)n 



(3.8) 

(3.9) 
(3.10) 



where {z; q)n and [z] q)oo are defined as 



1, 



if n = 0, 



(1 -z){l- zq){l - zq^) ■ ■ ■ (1 - ;zg"-^), if n > 



9)oo=n(i-^^")- 



(3.11) 
(3.12) 



n=0 



If we introduce the basic hypergeometric function [|T6 | 



2Vl 



^ (ai; q)n{a2] q)r 



n=0 



(6; q)niq;q)r 



(3.13) 
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( p.9| ) is rewritten as 



f{y) = y 



c 



(ay; q) 



2 

oo 



-2(Pl 



(3.14) 



This is the exphcit expression of the generating function 6(A,^). 

Next we are interested in the asymptotic behavior of Zl{C,) when L —>■ oo. This is 
determined by the singularity closest to the origin of the generating function B(A,^) as 
a function of A. We do not demonstrate this fact here; for a thorough discussion see 
Appendix A of ||12|. From the expression (|3.71 ), we see that there are two possible sources 

of singularities in 0(A,.^): a singularity of the function x (^if) (^^^ differentiation in the 

numerator does not change the position of any singularities) and a zero of the denominator. 
As we will see, the singularity which is closest to the origin is always of the latter type, for 
the region of the parameter space of interest here. 

First we consider the singularities of x(^j^^- From (|3.1CI|) it follows that y{l/u) is 
analytic at m = 0, that the singularity closest to the origin is the square root singularity 
at M = 1/4, and that \y{l/u)\ < 1 for \u\ < 1/4. Moreover, from ( p.9|) , the singularity of 
f{y) closest to the origin is a simple pole aX y = 1/a. Thus if Ao(0 denotes the singularity 

of X i^'ti) closest to the origin, then (i) Ao(0 = 1/(40 is a square root singularity if a < 1 
and (ii) Ao(0 = a/[(ct + 1)^^] is a simple pole if a > 1. In what follows we restrict our 
detailed analysis to the case a < 1. (See remarks following (|2.5| ); the analysis of a > 1 is 
similar.) 



In order to discuss the zeros of the denominator of (|3.71 ) we make the following assump- 
tion: 



If < g < 1 then the function / satisfies f\y) > for < y < 1. 



(3.15) 



This assumption is key for our analysis; we have considerable evidence for its truth, which 
we will discuss shortly, but at the moment no proof. From (|3.8|), ( p.lO| ), and ( p.l5| ) it 
follows that, when ^ is positive, the function x(^) increases monotonically from to /(I) 
as A increases from to \oiO- Thus if ^ satisfies < ^ < .^max, where ^max = /(I); the 
equation 



X 



1 

a! 



(3.16) 



has a unique solution A(0 in the interval [0, Ao(0]) and the denominator of ( p.7|) vanishes 
at this point. Since A(0 is smaller than Ao(0; this is the singularity of the function B(A, 
which is closest to the origin when < ^ < ^max- This singularity is a simple pole, so that 
Zl{^) behaves asymptotically (L oo) as 



Zl{^) ~ const. 



A(OJ 



(3.17) 
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Thus from ( plQ]) and (|2Tl| ), 



J(0= lim Ji(0 = (l-g)eA(0, (3.18) 



L— >oo 



These expressions can be simphfied somewhat because the function y{x) is exphcitly 
invertible. If we denote the inverse function of f{y) for < ?/ < 1 by g{C)y 

^ = fiy)^y = r\o = 9iO for o<e<emax, o<y<i, (3.20) 



then the equation ( 3.16|) can be rewritten as 



Combining the above results, we find 



^(OHl-.)^^^, (3.22) 



P(0 - 1 



{I + 9(0)9(0 



(3.23) 



From (|3.15|) and ( p.20| ) it follows that g ) increases monotonically from to 1 as increases 
from to ^max and, since g'(0 = ^/ f'(9(0)y ^^at < 5''(^) < oo in this region. Then ( p.22|) 



and (|3.23|) imply that J and p are analytic functions of ^ throughout the range < ^ < 
and that as ^ increases from to ^max, J(0 increases monotonically from to (1 — g)/4 and 
p(0 increases from to 1/2. We expect also, and will assume in what follows, that p(0 
is a monotonic function throughout the range < ^ < ^max. It follows that we may invert 
the function p(0 as Op) thus obtain an analytic function J(p) = J(Op)) defined for 
< p < 1/2. This is the indication that there is no mixed/disordered phase transition in 
the GCE. 

Note that for this argument it is important that f'(y) be strictly positive not only for 
< y < 1 but also for ?/ = 1 (and hence g'(0 finite not only for < ,^ < ^max but also for 
^ = ^max); this guarantees the convergence of p(0 in ( |3.23|) to 1/2 as ,^ — i^max. 

To obtain J(p) we must first invert the function f(y) to obtain g(0, then invert the 
function p(0 to obtain the fugacity as a function of the density. It is not possible to carry 
out these inversions explicitly in general, and we must use numerical methods to obtain 
J(p). But there are special cases where one or both of the inversions can be done explicitly. 
Let us study these cases first. 



Case 1. q = 0. For this case, the formula for f(y) in ( p.9|) is greatly simplified since 
(z; 0)n = (z; 0)oo = 1 — z and the only n = term in the infinite series of (O) remains. 
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We find 



fiy) 
f\y) 



y 



- ay 

ayf 



(3.24) 
(3.25) 



Notice that the assumption (|3.15|) holds in this case, since a > — 1. The function f{y) is 
invertible, and 



9{0 



1 + 2ae - yi + 4ae 



The exphcit formulae for J(^) and p(^) are then given by 

2a\ 



1 + a2 + 2a(l + a)2^ - (1 - a^)^! + 4a^' 

a(l + a)e[(l + a)Vl + 4ae - (1 - a)] 



(3.26) 

(3.27) 
(3.28) 



Vl + 4a^[l + a2 + 2a(l + a)^^ - (1 - a2)Vl +4a^] ' 

These formulae are further simplified when a = (as pointed out in the introduction, 
setting a = q = reduces the AHR model to the two-species ASEP studied in P]). Now 

f{y) = y, giO = and 



[1 + 0''' 



From the easy inversion of ( p. 301) we obtain the known current density relation 

J = p(l-p). 



(3.29) 
(3.30) 

(3.31) 



Case 2. a = —1. This case corresponds to taking the limit g — > 1 with a fixed, and is 
thus strictly speaking outside our parameter region < g < 1, but it is nevertheless of 
interest. Moreover, it will play the role of zeroth-order approximation in the analysis of 
the next section; for this reason, the functions for this special case will be distinguished by 
the superscript (0) in the following. Using the formula for 20i ||16|| , 



2(Pl 



0-1, fl2 
b 



aia2 



(6/ai;g)oo(&/o2, ;«?)c 



(&;g)c 



' b 



we find 



y 



d_ 

dy 



f^'\y) 



yr 

- y 



(1 + ^/) 



3' 



(3.32) 

(3.33) 
(3.34) 
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Assumption ( |3.15| ) is not strictly satisfied, since df^'^\y)/dy = for ?/ = 1, but f{y) is 



strictly monotone for < y < 1 so that the derivation of the asymptotic form (|3.17|) is 
still valid. Note that ( |3.33| ) and (|3.34|) are independent of q and hence may be obtained by 



setting a = — 1 in ( p.24| ) and ( p.25|) . Hence the explicit formulae for g{^), J{^) are obtained 



by simply taking a = — 1 in ( p.26| ) and ( p.27| ): 



,<«.(,) . i^ii^^. (3,35) 

j(°)(0=e (3.36) 

As for p(^), we have p^°^(0 = for < ^ < ^max = 1/4, while p(°)(l/4) is not well specified. 
Thus we cannot achieve nonzero density values by any specification of the fugacity: the 
q ^ 1 limit is a very singular limit for the density-fugacity relation. 

We now turn to the question of the validity of the key assumption ( p.l5|) . We believe. 



although we have not proved, that (|3.15|) holds for all a,q satisfying — 1 < a < 1 and 
< g < 1. This belief is based on several pieces of evidence. We have evaluated the 
functions f{y) and f'{y) numerically for various values of a, q, and y and have always 
found that /'(|/) > for < ?/ < 1 (such evaluations are easily carried out to arbitrarily 
high accuracy from the formula ( p.9|) ). We do observe, however, that for q and y close 
to 1 the value of f'{y) is very small; for instance, /'(I) is of order 10~^^ when q = 0.9 
and a = 1. This phenomenon will play a key role in the analysis of section 4. Note also 
that it follows from the considerations of the special cases above that ( |3.15| ) holds when 
g = 0, and that when a = — 1 (recall that this corresponds to a limit q ^ 1), f'{y) > 
when < y < 1; this suggests that f'{y) > except at the corner q = y = 1 of the 
region of interest. Moreover, we will show in the next section that /'(I) > when a = 1 
(i.e., q = —a). This is a weaker statement than (|3.15|) but, since f{y) is expected to be a 



monotonically decreasing function, is a strong analytical indication that ( |3.15| ) holds when 
q = —a. 



4 Two Distinct Regimes without Non-analyticity 

In the previous section we argued that, even in infinite volume, the current is an analytic 
function of the density throughout the mixed/disordered phase. On the other hand, it is 
clear from inspection of the J-p diagram in the thermodynamic limit that, for large values 
of q, the derivative of J(p) does change rapidly — in fact, extraordinarily rapidly — as p 
passes through some critical value pc- For example, when q = —a = 0.9, there is a change 
of order one in J'{p) arising from a change of order 10~^^ in p (see Fig. 2). In this section, 
we explain how this behavior arises: it is associated with the very small values of /'(I) 
occurring when q is near 1, which lead to rather peculiar behavior in the density-fugacity 
relation ( p.23| ). 

In Fig. 3, we graph the functions J(^) and p(^) of (|3.22|) , ( p.23| ) for a = 1 and several 



values of q. Note that, although the function J(^) remains quite smooth as q increases, and 
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indeed even for q = 1, the function p(^) develops an apparently sharp corner or kink, even 
for values of q about 0.75; this kink corresponds to the kink in the J-p curve at Pc- Here we 
give an approximate analysis of this behavior under the assumption that q is reasonably 
near to 1 and that /'(I) is very small. 

We first analyze the region p < pc. Using the formula (see |T6[) 



2Vl 



ai, 02 



(02; q) {aiz; q) 
{b; q)ooiz;q)oo 



-2HJ1 



6/02, z 



(4.1) 



we obtain an alternate expression for f{y): 

I 

y 



(4.2) 



Since we are interested in values of q close to 1 we set g = 1 — e, where e(> 0) is assumed 
to be small, and expand in e, denoting quantities accurate to order /c by a corresponding 
superscript (e.g. /'•'^•*). The 0**^ order was analyzed in the previous section; here we will 
carry out the expansion to first order. By a straight-forward computation we obtain 



f'\y) 



[i + yf 



2y 



f^'\y) 



V^\y) 



a 



(4.3) 



As y increases from zero to one, f^^\y) increases monotonically from zero to .^1 
1/4 + e/iSa). The function f^^\y) is invertible, with inverse 



(1) 

max 



iO)( 



1-2X(0- v/l-4X(e) 
2X(0 



(4.4) 



where 



x{0 



a 



a 



(4.5) 



for < ^ < ^mlx- Then the density-fugacity relation (|3.23| ) becomes 



p(0 



so that to first order in e the density is 



x{0 



a 



0(e' 



(4.6) 



(4.7) 



Comparing these results with the curve for g = 0.75 in Fig. 3, we realize that ( |4.6| ) gives 
the behavior of the density-fugacity relation up to the kink; since this region includes ^ = 
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(or p = 0), it corresponds to the disordered regime. As the fugacity ^ increases from zero 
to .^max, the density p*^"*^^ increases hnearly from zero to e/(4a), so that this value may be 
identified, up to first order, with the density pc- 

Pi" - (^•«> 

Note that under this analysis there is no value of ^ which gives rise to a density p in the 
interval (pc, |]. This is because our approximation of /(?/) by f^^\y) lacks an important 
property of the original f{y): the derivative of f^^^ becomes zero at ?/ = 1, whereas our 
basic assumption ( |3.15| ) is that /'(I) > 0. In fact, it appears that /'(I) vanishes to all 
orders in the perturbation expansion in e; this is shown below (see (|4.24 )) in the case 
q = —a. 

We now analyze the region Pc < p ^ \ under the assumption that /'(I) is positive but 
very small. We will also assume (as is supported by numerical evaluations) that /"(I) is 
negative and of order unity. Under these assumptions, and recalling that /(I) = ^max, we 
may approximate /(?/) near y = 1 as 

f{y) ^ r{y) = U.. - B{{1 + 6-y)'- 6'), (4.9) 

where B = —/"(I) and 6 = — /'(!)//"(!); that is, f{y) will have a quadratic maximum at 
approximately 1 + 6, where 6 is very small. The function f*{y) has inverse 



g*{^) = 1 + 6- ^^-^^^^^ + 6^ = 1 + 6-6./^, (4.10) 
where for ^ ~ ^max, it is convenient to introduce the scaled variable rj defined by 

emax-e = ^<5'- (4.11) 

The density-fugacity relation ( |3.23| ) then becomes p(^) ^ P*{Oj where 



^max + 1 



4B + T 
so that as r] goes infinity, p*(0 approaches 

. 1 



(4.12) 



^^ = 2 



1 



45 



(4.13) 



Thus here the critical density pc arises as the minimum value of the density which can 
be described by (|4.12| ). Clearly (|4.12| ) describes the density-fugacity relation in the mixed 



regime, that is, for p > pc- When q = 0.75, for example, this corresponds to the nearly 
vertical line shown in Fig. 3. The value pc of (|4.13 ) is an extremely accurate value for the 



"critical density" between the mixed and disordered regimes when q is close to one; this is 
the value used in plotting Fig. 2 and thus in the case q = —a = 0.9 is accurate to about 

10-24. 
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Note that ( |4.3| ) implies that to first order in e, B = 1/16 + e/(16a); recalhng that 
d'lx = 1/4 + e/(8 a), we see that ( [4 .131 ) reduces to ( |4.8| ) in first order. In the critical 
value of p was estimated, on the basis of extrapolation of finite size data, to be (in our 
notation) (1 — g)/[4a — 2(1 — q)]; this value agrees with ( [4.13| ) to first order in e. 

The key point for the above analysis, as we have emphasized, is that /'(I) is "very 
small" over a range of q values close, but not necessarily very close, to 1. We can throw 
some additional light on this phenomenon because, in the special case a = 1, it is possible 
to compute /'(I) exactly. When a = 1, i.e., when a = —q, the expression of f{y) in ( |4.2| ) 
can be further simplified as 



fiv) 



n=0 



yq 



(4.14) 



so that 



/'(I) 



4 



„ (1 + g")2 2^ iJ^qn 



n{—qY 



+ q'' 



(4.15) 



_n.=0 n=l 

The sums in ( [4. 151 ) can be evaluated explicitly, using the formulae [|I6| 

1 4 



^ (1 

n=l ^ 



-qr 



+ q 



n\2 



n=l 



+ q'- 



{q; q)oc 
[-q; q)c 



1 



(4.16) 



This leads to 



/'(I) 



1-g 



{q\q)oc 
[-q; q)c 



1 4 



(4.17) 



This expression is valid for any value of q satisfying < g < 1, and is clearly positive. The 
behavior of /'(I) function of q is shown in Fig. 4. 

Now we consider the asymptotic behavior of /'(I) as e = 1 — g becomes small. Setting 
q = e~'^ *, we find 



{q\q)oc 
[-q; q)c 



n 

n=l 



1-q' 



^9o(0,ivrt). 



(4.18) 



Here 'do{x, r) is an elliptic theta function defined by 

oo 

^o(a;,r) = 



n ^inrn^ +2mxn 



(4.19) 



Applying Jacobi's imaginary transformation, 

1 



(4.20) 
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with 

oo 
n=— oo 

we obtain 



-g;g)oo V^rt 



^ g-(n-l/2)Vt. (4_22) 



n=— oo 



(These formula involving theta functions are easily obtained from, for example, Chapter 21 
of 1^.) Now it is easy to read off the asymptotic behavior as q approaches 1. The leading 
contributions come from the n = and n = 1 terms in ( |4.22|) . Since t ~ e/vr^ + e^/(27r^), 
we find 

- ^ e-^/^*^2./^e^e--^/^^ (4.23) 



(-9;?)oo V^rt V e 

Finally we obtain the asymptotic expression for /'(I) 



/'(I) ^ ^e^^ e-"'/^ (4.24) 



Using ( 3.23 ), it is also possible to obtain the asymptotic expression for 5: 



, /'(I) 2V .2 

5= ~ ^e-e"" /^ (4.25) 

2r) e 

5 Density Fluctuations and the Canonical Ensemble 

In the previous two sections, we have studied the current in the thermodynamic limit in 
the framework of a GCE. Our analysis strongly suggests that there is no phase transition 
between the mixed and disordered regimes in the GCE. In this section we argue that in 
the thermodynamic limit, the current in the GCE and the CE are the same, that is, that 
the current J{p) discussed in Section 3 and 4 is the same as the infinite volume limit of the 
current in the CE ensemble with density p. This implies that the lack of a phase transition 
in the GCE discussed above also holds for the CE. The argument also sheds some light on 
the structure of the infinite volume state of the model in the mixed/disordered regime. 

In finite volume, the GCE is a superposition of canonical ensembles: contributions come 
from all values of the densities, with weights which are determined by the fugacity ^. The 
total mean density of particles (both + and — ) is 2pi(^), where Pl{C) is given by formula 
( p. lip . The fluctuations in the total density can be expressed similarly (ideally we should 



discuss fluctuations in the density of each species, but that seems more difficult). If we 
write = + N_ and introduce the quantity 



m) = \{{N')L-{N)l) = \{il^\nZ, 



(0, (5.1) 
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then the fluctuations in the total density are simply given by Fl{^)/L. Now if we take the 
thermodynamic limit in the GCE for some fixed ^, then the mean total density 2pl{^) is 
expected to have a well defined limiting value, 2p(^) = limi^oo '2pl{,C)- If we can show that 
the fluctuations in the total density go to zero in the thermodynamic limit, then we expect 
that in this limit only a single total density 2p((^) survives, and if we ignore the possibility 
of canceling fluctuations in the + and — densities then these densities will also be unique 
in the limit. Thus, in the thermodynamic limit, the GCE with fugacity ^ is the same as 
a CE with (+ or — ) density This corresponds to the usual GCE/CE equivalence for 
equilibrium systems away from a phase transition. For our model, of course, condensation 
of particles in the mixed regime suggests a possible phase transition, and we need an 
independent argument. It turns out that we can exactly calculate F{^) = lim^^^oo -^^(0 
from (lOSi ): 



- C(i - gmgiOg'iO - (i - g^giOgiO ) ■ (5.2) 



Numerical values of this quantity can be very large. For instance, for q = —a = 0.9, we 
find that F{C,) is of order 10™ in the region p > p^. However, the important point is that 
Fl{^) has a finite value in the thermodynamic limit, so that the density fluctuations F^/L 
in the finite system vanish in the infinite system. As argued above, this suggests that our 
analysis in the previous two sections, of the behavior of J{p) in the GCE, is valid for the 
CE as well. 

It is also possible to estimate the length scale R up to which particle condensation will 
be observed in Monte Carlo simulations. Let us write {N"^) = Ylii<i jKi^^i^j) (|5TT|), 
where = 1 if there is a particle at site i and aj = otherwise. If we assume that the 
particle-particle correlation function decays exponentially with some characteristic length 
R, {aiaj) ~ e-l^--'!/^, then we see that Fl ~ R. Thus the correlation length can become 
huge; for instance, for the above case of g = —a = 0.9, R would be of the order 10''°. 
Therefore it is certainly not possible to observe the breakdown of the particle condensation 
with Monte Carlo simulations. 

The discussion of this section also throws light on our assumption of section 3, that 
p{^) is a strictly increasing function. Since 2^dpL/d^ = is a fluctuation, Opl/O^ must 
be positive, so that dp/dC, > 0. Vanishing of dp/dC, could occur only with subnormal 
fluctuations in (lim^^^oo -^l = 0), a scenario which is certainly unlikely and has been 



ruled out in certain equilibrium situations ||26| 



6 Concluding Remarks 



In this article we have studied an exclusion process on a ring, originally introduced and 
investigated by Arndt et al []I1|,|TH|. Using the theory of g-orthogonal polynomials, we have 



obtained the exact expression of the current in the thermodynamic limit in the framework 
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of a grand canonical ensemble (GCE). Contrary to what is suggested by Monte Carlo 
simulations and finite volume calculations, [0,11^], we find no phase transition between 
the mixed and the disordered regime; specifically, the infinite volume current as a function 
of the density is analytic through the mixed/disordered phase. However, the derivative of 
the current can change rapidly over a very small (but finite) interval of p, especially if q is 
close to 1. 

Our analysis of the infinite volume current depends on a key assumption, inequality 
( p.l5| ). Although this is unproved, we have strong evidence of its truth: it holds in the 
special case g = (and essentially also the special case a = —1); it is supported by the 
analysis of section 4 for the case q = —a; and we have checked its validity for various 
parameter values by highly accurate numerical computations. However, the rigorous proof 
of ( p. 151 ) remains an outstanding problem. A second assumption, the positivity of p'{C,), is 
physically natural and is also supported by numerical calculation. 

Most of our analysis was carried out in a GCE. The relation between the canonical 
and grand canonical ensembles in equilibrium statistical mechanics has been well studied, 
but we know of no general results for non-equilibrium situations. Here we have simply 
asked how far the results for our GCE (which itself represents a particular and perhaps 
non-canonical choice) are also valid for the CE. We have given a heuristic argument: exact 
calculations show that the density fluctuations in the GCE vanish in the thermodynamic 
limit, which suggests that, in this limit, the GCE for a given fugacity agrees with the CE 
for a certain corresponding density. 

It is fascinating that in the relatively simple AHR model one already has a quite subtle 
transition between the disordered and the mixed regime. As discussed in section 5, com- 
puter simulations would need to be carried out up to lattice sizes of the order 10^° (when 
q = —a = 0.9) to see the breakdown of the particle condensation, so that this is certainly a 
phenomenon which one would not expect from the accessible finite size Monte Carlo data. 
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Figure Captions 



Fig. 1: Space time diagrams of the AHR model from Monte Carlo simulations. The 
horizontal axis represents the site number j whereas the vertical axis represents time. The 
existence of a positive (resp. negative) particle is represented as a black (resp. gray) point. 
For all three figures, we set L = 100, = A^„ = 30, a = 0.5. The left-most figure 
corresponds to the pure phase {pr = 0.9, pl = !)• The figure in the middle corresponds to 
the mixed regime [pfi = 1,pl = 0.9). The right-most figure corresponds to the disordered 
regime {pp = 1,pl = 0). 

Fig. 2 : The J-p diagram for the case a = l,q = 0.9. The horizontal and vertical axes of 
the inset are (p — pc) x 10^^ and (<//(! — q) — 0.25) x 10^^ respectively. The function J(p) 
is analytic, although its derivative changes very rapidly as p passes Pc- 

Fig. 3 : The functions J(^) and p(^) when a = 1. For g = 0, we have J(^) = + 
(PI) and p(0 = + (Pi- For g = 1, we have J(0 = ^ (P|) and p(0 = (0 < 
^ < 1/4). The function J(^) for < g < 1 interpolates these two limiting cases smoothly. 
On the other hand, as q increases from zero, there appear two distinct regions for p{C,)- In 
particular, the g — > 1 limit appears to be singular. 

Fig. 4 : The behavior of /'(I) as a function of g. We see that /'(I) becomes extremely 
small as g approaches one. 
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